Draft version June 20, 2008 

Preprint typeset using I4TgX style emulateapj v. 1 1/26/04 



THE EVOLUTION OF THEM BH -er RELATION INFERRED FROM THE AGE DISTRIBUTION OF 
LOCAL EARLY-TYPE GALAXIES AND AGN EVOLUTION 

Francesco Shankar 1 , Mariangela Bernardi 2 , and Zoltan Haiman 3 

Draft version June 20, 2008 



00 

o 

O 

(N 

a 

>—> 
O 

6 



ABSTRACT 

We utilize the local velocity dispersion function (VDF) of spheroids, together with their inferred age- 
distributions, to predict the VDF at higher redshifts (0 < z < 6), under the assumption that (i) most of the 
stars in each nearby spheroid formed in a single episode, and (ii) the velocity dispersion a remained nearly 
constant afterward. We assume further that a supermassive black hole (BH) forms concurrently with the stars, 
and within ± 1 Gyr of the formation of the potential well of the spheroid, and that the relation between the mass 
of the BH and host velocity dispersion maintains the form Mbh oc a 13 with /3«4, but with the normalization 
allowed to evolve with redshift as oc (1 + z) a - We compute the BH mass function associated with the VDF at 
each redshift, and compare the accumulated total BH mass density with that inferred from the integrated quasar 
luminosity function (LF; the so-called Soltan argument). This comparison is insensitive to the assumed duty 
cycle or Eddington ratio of quasar activity, and we find that the match between the two BH mass densities fa- 
vors a relatively mild redshift evolution, with a ~ 0.26, with a positive evolution as strong as a > 1 .3 excluded 
at the 99% confidence level. A direct match between the characteristic BH mass in the VDF-based and quasar 
LF-based BH mass functions also yields a mean Eddington ratio of A ~ 0.5 — 1 that is roughly constant within 
< z < 3. A strong positive evolution in the Mbh-c relation is still allowed by the data if galaxies increase, on 
average, their velocity dispersions since the moment of formation, due to dissipative processes. If we assume 
that the mean velocity dispersion of the host galaxies evolves as <r(z) = c(0) x (1 + z) ',we find a lower limit 
of 7 > 0.23 for a > 1.5. The latter estimate represents an interesting constraint for galaxy evolution models 
and can be tested through hydro simulations. This dissipative model, however, also implies a decreasing A at 
higher z, at variance with several independent studies. 

Subject headings: : black hole physics - galaxies: evolution - galaxies: active 
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1. INTRODUCTION 

It has now been assessed that most, if not all, local galaxies 
have a supermassive black hole (BH) at their center, the mass 
of which is tightly correlated with the velocity dispersion a 
and other bulk properties of the host galaxy (e.g., Ferrarese & 
Merritt 2000; Gebhardt et al. 2000). However, the sample of 
local galaxies for which the BH sphere of influence has been 
resolved amounts to only ~ 30. It is not clear how represen- 
tative this small sample is of the whole BH population, and 
whether the correlations seen in the sample already held in 
the past. 

Peng et al. (2006) have collected a sample of 31 lensed 
and 1 8 non-lensed Active Galactic Nuclei (AGN) at redshifts 
z > 1.7. They measured rest-frame /?-band luminosities from 
//-band fluxes and BH masses by applying virial relations 
based on emission line widths. They found that the BH-to- 
host galaxy luminosity at z ~ 2 is about the same as the one at 
z ~ 0. Therefore, once the observed rest-frame luminosity is 
dimmed through passive evolution to z ~ 0, at fixed BH mass 
the ratio BH-to-host luminosity grows significantly, and the 
resulting BH-luminosity normalization is several times higher 
than the local one. Similar results were derived by McLure et 
al. (2006), who measured the BH-to-host galaxy mass ratio in 
a sample of radio-loud AGNs in the redshift range < z < 2 



finding M B h/M star oc (1 



Shields et al. (2006) found 



1 Department of Astronomy, The Ohio State University, Columbus, OH 
43210; shankar@astronomy.ohio-state.edu 

2 Department of Physics and Astronomy, University of Pennsylvania, 209 
South 33rd St, Philadelphia, PA 19104; bernardm@physics.upenn.edu 

3 Department of Astronomy, Columbia University, 550 West 120th Street, 
New York, NY 10027; zoltan@astro.columbia.edu 



that the CO emission lines in a sample of z > 3 quasars is 
very narrow, suggesting bulge mass about an order of mag- 
nitude lower than measured in the local universe, at fixed 
BH mass (see also Coppin et al. 2008). Treu et al. (2007) 
found that the BH masses in a sample of 20 Seyferts galax- 
ies at z = 0.36 are offset by an amount of AlogMeH^ 0.5 
at fixed velocity dispersion, which implies an evolution of 
Mbh/Mstaro^ (1 + z) 1,5±L0 , consistent with that derived by 
the previous works. 

On the other hand, Lauer et al. (2007) have discussed sev- 
eral possible biases which may seriously affect these findings. 
At high redshifts a sample will be biased toward the most lu- 
minous AGNs and more massive BHs. Given the observed 
scatter in the local relations, especially significant in the Mbh- 
host luminosity relation, these massive BHs will be preferen- 
tially associated with the less massive, but more numerous 
galaxies, yielding a false sign of evolution. When the cumu- 
lative mass density of AGNs is taken into account, several 
authors (e.g., Haiman et al. 2004; Marconi et al. 2004; Sil- 
verman et al. 2007; Shankar et al. 2008a, hereafter SWM) 
have shown that once rescaled by a simple constant, it pro- 
vides a good match to the cosmological star formation rate 
density. De Zotti et al. (2006) and SWM have shown that 
the galaxy stellar mass function at z ~ 2, mostly composed 
of massive early-type galaxies (e.g., Drory et al. 2005), con- 
verted into a BH mass density assuming a Mbh/A^star ratio 
3-5 times higher than in the local universe, would imply a BH 
mass density already close, if not higher, than that inferred 
in the local universe, leaving no room for further accretion at 
z < 2, where, in fact, a significant fraction of the total AGN 
energy output is produced. Recently, Ho et al. (2008) com- 
piled a sample of 154 nearby (z < 0. 1) active galaxies showing 
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substantial ongoing BH growth in the most actively accreting 
AGNs, where BH growth appears to be delayed with respect 
to the assembly of the host galaxy. 

In this paper, we propose a simple, yet robust, way to con- 
strain the degree of redshift evolution in the Mbh-c rela- 
tion, that is relatively insensitive to assumptions that relate 
the SMBH population to quasars. We combine the measured 
velocity dispersion function (VDF) of local spheroids with a 
postulated power-law redshift-dependence of the Mbh-c re- 
lation. By comparing the resulting total BH mass density at 
each redshift with the same quantity inferred from integra- 
tion of the active galactic nuclei (AGN) luminosity function 
(see Soltan 1982), we find the degree of evolution required 
in the Mbh-c relation to match these two independent es- 
timates. This approach yields results based on the "bulk" 
of the active BHs at all redshifts, and is therefore relatively 
insensitive to possible biases which may affect studies per- 
formed on small samples of high-redshift luminous quasars 
(e.g., Lauer et al. 2007). After describing the sample used in 
our computations in § we proceed to derive our main re- 
sults in § |3] These results are discussed further in § [5] where 
we also offer our conclusions. Throughout this paper we 
use the cosmological parameters fl m — 0.30, 57a = 0.70, and 
h = i/o/lOOkms^'Mpc -1 = 0.7, consistent with the three- 
(Spergel et al. 2007) and five-year (Dunkley et al. 2008) data 
from the Wilkinson Microwave Anisotropy Probe (WMAP). 

2. DATA 

We have used the sample of early-type galaxies obtained 
by Bernardi et al. (2006). The sample, extracted from the 
Sloan Digital Sky Survey (York et al. 2000), contains over 
40,000 early-type galaxies, selected for having an apparent 
magnitude 14.5 <M, < 17.75, extending over a redshift range 
0.013 < z < 0.25, which corresponds to a maximum lookback 
time of 3 Gyr. The ages of galaxies are computed in two dif- 
ferent ways, discussed in detail by Jimenez et al. (2007), from 
(i) single stellar population spectral fitting, using the MOPED 
algorithm (Heavens et al. 2000) to determine the full star- 
formation history of the galaxies, and (ii) using the published 
ages by Bernardi et al. (2006) which were obtained by fitting 
the Thomas et al. (2005) a-enhanced models to the Lick index 
absorption features measured from stacked spectra of galaxies 
with similar properties. The age-distributions at fixed veloc- 
ity dispersion a are generally broad, but tend to be narrower 
and centered on older ages for higher values of a. Such effects 
are more marked for age distributions inferred from MOPED 
(see Figure 1 in Haiman et al. 2007). We will compare results 
obtained by adopting either the MOPED or the Lick-index age 
distribution in § [3] 

The analysis presented in Bernardi et al. (2006) and 
Haiman et al. (2007) probe velocity dispersions within 2.05 < 
log(cr/kms _1 ) < 2.45. Here we extend such analysis in- 
cluding the age distributions of galaxies with velocity dis- 
persion 2.45 < log(CT/kms _I ) < 2.55. We find that galax- 
ies within this last bin are even older than the oldest galax- 
ies probed by Haiman et al. (2007), confirming and extend- 
ing the general trend of increasing age for larger a. Instead 
of considering one single bin with mean velocity dispersion 
log(cr/kms _1 ) = 2.5, we have treated the bin as two distinct 
bins with log(cr/kms _1 ) = 2.45 and 2.55 which we have as- 
sumed share the same age distributions as the total bin. We 
have also included an additional bin with log(cr/kms _1 ) = 
2.60, which we have again assumed to have an age distribu- 



tion equal 4 to the one with log(cr/kms _1 ) = 2.45 and 2.55 
(a direct estimate of the ages for these galaxies with the tech- 
niques discussed above is highly limited by the low Signal- 
to-Noise of the spectra). As will be shown in §[3] (Figure [T), 
this binning in log a enables us to better probe the statistical 
evolution of the VDF even at large velocity dispersions, and it 
has a negligible effect in the resulting cumulative black hole 
mass density and on our general results. 

3. RESULTS 

We first estimate the VDF as a function of redshift z. At any 
z, the VDF in a given bin of velocity dispersion er, is given by 
all the galaxies which have formed prior to z. Therefore to 
compute $(ct/,z) we subtract from the local census of galax- 
ies with velocity dispersion 07 those galaxies that have an age 
t lower then the lookback time 77 (z), 



i-X>fa(z)ki) 



T<Tj 



x $(07) . 



(1) 



Note that p(t/(z)|<7,) refers to the fraction of galaxies with ve- 
locity dispersion cr, which have an age of t,-(z)± 1 Gyr. There- 
fore <&(<T(,z) includes in the ct, bin all galaxies whose ages are 
within ± 1 Gyr of 77 (z). The VDF at z = is taken from Sheth 
et al. (2003) and includes the contribution of bulges of spirals. 
We therefore assume that bulges of spirals and local spheroids 
within the same bin of velocity dispersion share similar age 
distributions. However, as discussed in §[5]below, our results 
would still hold even if the contribution from spirals were ne- 
glected. The statistical uncertainties associated with $(cr,z) 
are computed from equation (Q~|) through error propagation in- 
cluding uncertainties in p(Tj(z)\<r), given by Haiman et al. 
(2007; see their Figure 1) and 3>(<r), given by Sheth et al. 
(2003). 

Figure Q] shows the VDF obtained from equation (Q~|). The 
different curves in both panels show Q(o-,z) at different red- 
shifts, as labeled. At fixed redshift, the symbols indicate the 
position of the mean in the bin of log a considered, for which 
reliable age distributions »(t) (z)|0';) have been computed. In 
the left panel of Figure [T] the pj, distributions have been de- 
rived from the Lick-indices method, while the right panel 
shows the results with the MOPED-based py, distributions. 
In our analysis below, we will adopt the Q(a,z) implied by 
the py derived from Lick indices. However, we will discuss 
the consequences of the alternate choice on our results in §[5j 
In both cases, we reproduce the conclusion of previous work 
(e.g., Trager et al. 2000; Thomas et al. 2005; Bernardi et 
al. 2006; Jimenez et al. 2007; Haiman et al. 2007) - that is, 
we find strong evidence for downsizing: on average, galaxies 
with larger velocity dispersion are formed earlier. This be- 
havior is expected from basic galaxy formation theory: high- 
redshift galaxies form in a denser universe and therefore pref- 
erentially form out of baryonic clumps collapsed in denser, 
gas rich environments which in turn, induce more dissipation, 
more compact remnants and higher velocity dispersions. At 
fixed velocity dispersion, the MOPED ages are higher than 
infeiTed from Lick indices, producing a less pronounced evo- 
lution in the VDF at < z < 3. Theoretical models in which 

4 a more appropriate choice would be to assign older ages to the galaxies 
with extreme velocity dispersions, given the general trend of older ages for 
higher a, however this would pose even stronger evidence for the downsizing 
discussed below further strengthening our general conclusions. 
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FIG. 1 . — Velocity dispersion function (VDF) at different redshifts, as labeled, obtained by combining the local VDF with the age distribution of local galaxies 
(equation [I]). Left panel: VDF obtained using the ages computed from Lick indexes; Right panel: VDF obtained using the ages computed from the MOPED 
algorithm (see text for details). 



the galaxy velocity dispersion is linked with the virial veloc- 
ity of the host halo (e.g., Ferrarese 2002) predict similar trends 
for the VDF as a function of time (Cirasuolo et al. 2005; see 
also Loeb & Peebles 2003). 

The BH mass function implied by the VDF at any time is 
given by converting <P(cr,z) to a BH mass function through 
the Mbh-c relation and a convolution with a Gaussian with 
intrinsic scatter of 0.22 dex. We assume here that logMeH at 
fixed logcr is given by a Gaussian distribution, with a mean of 



log 



M 



BH 



MQr 8.21 + 3.831o gl ^-^j + «log(l + z), 

(2) 

and a standard deviation of r\ = 0.22. This latter value rep- 
resents the intrinsic scatter as given by Tundo et al. (2007) 
and as recently confirmed by Shankar & Ferrarese (2008). By 
integrating the resulting BH mass function at all times, we de- 
rive the total BH mass density Pvdf(z), corresponding to BHs 
in the range of a and BH mass probed by our sample at each 
redshift. 

We then compare Pvdf(z) with the BH mass density ob- 
tained by direct integration of the AGN luminosity function 
<3?(L,z) from z = 6 up to redshift z. The latter quantity is given 
by 

1 — e r 6 r°° dt 
p. (> log L min ,z) = — T dz' §(L,z')L—d\ogL. 

(3) 

Here e represents the radiative efficiency, and for our numer- 
ical calculation, we adopt the bolometric AGN luminosity 
function $>(L,z) from SWM (using the LF from Hopkins et 
al. 2007 gives similar results as discussed below). At each 
redshift, we integrate equation © above the minimum lumi- 
nosity, corresponding to the minimum BH mass probed by 
the velocity dispersion distribution in our sample via equa- 
tion (fJJ, which is MBH.min ~ 10 7 M Q , in the reference model. 
The minimum luminosity is then taken to be the Eddington 
(1922) luminosity corresponding to this BH mass. As we dis- 
cuss in § |5] below, the exact choice of the lower limit in the 
integral of equation (f3]i does not alter our conclusions. 

The growth rate of an active black hole of mass Mbh is then 
A?bh =Mmi/t e f, where the e-folding time is (Salpeter 1964) 



fef =4X 10 7 



e(l-e) 



0.1 



A-'yr, 



(4) 



where A is the ratio of the luminosity eM^c 2 to the Edding- 
ton luminosity. Figure [2] compares the two independent es- 
timates of BH mass densities. The accreted mass density at 
each redshift obtained from <f>(cr,z) and the Mbh-c relation 
is shown with a solid curve. The solid squares show the red- 
shifts where the mass density was computed. The long-dashed 
curve represents the mass density inferred from integration of 
the SWM AGN bolometric luminosity function. Given that 
the ages of galaxies in the sample have a median associated 
uncertainty of ± 1 Gyr, at any time t(z) the BH mass density 
from AGNs to be compared to /Ovdf(z) is systematically un- 
certain by the mass accreted within t ± 1 Gyr, which we show 
as the gray area. 5 We choose a constant mean radiative effi- 
ciency of e = 0.080, which provides a good match to the BH 
mass density at z — (e.g. Haiman et al. 2004, SWM). It can 
be immediately inferred from the left panel, which assumes an 
unevolving Mbh-c relation, that /Ovdf(z) and p,(z) are con- 
sistent with each other within errors, and therefore a strong 
evolution with redshift in the Mbh-o" relation is not required. 
Very similar results are found if we adopt the bolometric lu- 
minosity function from Hopkins et al.(2007), shown as the 
dot-dashed curve in the same Figure. In this case, we use a 
slightly higher radiative efficiency of e = 0.104 to renormal- 
ize the total z = accreted mass density to the local value, 
due to the fact that the bolometric corrections used by Hop- 
kins et al. (2007) are about 30% higher then those adopted by 
SWM. Nevertheless, even in this case we find that p, (z) well 
matches Pvdf(z) at all times. 

Joint confidence levels on the two parameters e and a, in- 
ferred from a x 2 analysis are shown in Figure [3] The cross 
marks the best-fit model with e = 0.08 and a = 0.26 (corre- 
sponding to the minimum x^in ^2.6 for 8 degrees of free- 
dom), which is shown in the right panel of Figure [2] Once 
a constant radiative efficiency is fixed to match the z = lo- 
cal and accreted mass densities, it is evident that the available 
data favor a relatively mild redshift evolution of the Mbh-c 
relation with a < 0.3, while a strong evolution with a > 1.3 
is ruled at 99% confidence level. Likewise, negative evolu- 
tion with a <y — 1 is ruled out for any choice of e. We note 

5 Note that the t ± 1 Gyr uncertainty is for Pvdf(z)- However, in our calcu- 
lations, assigning the uncertainty to p. (z) or Pvdf(z) makes no difference. If 
the time of formation of the galaxies is uncertain by ± 1 Gyr, then statistically 
the pvDF(z)at the time 1 can be compared with the cumulative mass accreted 
at any time (±1 Gyr. 
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FIG. 2. — Comparison between the accreted mass density at each redshift obtained from $(cr,z), convolved with the local Mgn-o" relation (solid curve and solid 
squares) and the mass density inferred from integration of the Shankar et al. (2008a) AGN luminosity function {long dashed curve) and a radiative efficiency of 
e = 0.080; the grey area represents the uncertainty at each time t associated to the mass accreted within t ± 1 Gyr; the dot-dashed line is the predicted accreted 
mass using the Hopkins et al. (2007) AGN luminosity function and a radiative efficiency of e = 0.104. The left panel shows Pvdf(z) predicted from the VDF 
assuming that the Mbh-o" relation is independent of redshift, while the right panel shows the predictions for the best-fit Pvdf(z) when the normalization of the 
M-gw-o- relation evolves as oc (1 +z) 26 - 



that values of a > 0.8, in fact, yield the unphysical result that 
the absolute total BH mass density increases from z = to 
Z > 0.7, as shown in Figured The confidence contours shown 
in Figure 0may therefore overestimate the maximum allowed 
value of a. 

It is clear from equation (0 that the accreted BH mass 
density does not depend on the assumed duty cycle or Ed- 
dington ratio distribution A(Mbh,z), apart from a weak de- 
pendence on the latter through the lower limit of the inte- 
gration. The strongest dependencies are on the radiative ef- 
ficiency and on the bolometric corrections (see also Figure 9 
in SWM). On the other hand, the Eddington ratio distribution 
and its evolution with redshift, can be constrained by compar- 
ing the AGN-based and VDF(z) -based differential BH mass 
functions (rather than comparing only the integrated quanti- 
ties). The $(cr,z) convolved with the Mbh-c relation (equa- 
tion[2ji in fact, predicts the shape of the BH mass function for 
Msh ii A^BH.min- On the other hand, as extensively discussed 
in the literature (see SWM, and related work by, e.g. Cavaliere 
et al. 1982; Small & Blandford 1992; Salucci et al. 1999; Yu 
& Tremaine 2002; Marconi et al. 2004; Shankar et al. 2004), 
if a mean Eddington ratio A = L/Ledd is assumed for the ac- 
tive BHs, then through a continuity equation and an assumed 
initial condition, the AGN luminosity function can be directly 
mapped into a BH mass function at all times. The "break" 
in the predicted BH mass function will then approximately 
reflect the break L* (z) in the observed AGN luminosity func- 
tion, i.e., Mbh*(z) ocL*/A, where A is the mean Eddington ra- 
tio. Following Soltan (1982) and Salucci et al. (1999), SWM 
(see also, e.g., Yu & Tremaine 2002) showed that constraints 
on the mean radiative efficiency and Eddington ratio of BHs 
can be gained by comparing the directly measured and the ac- 
creted BH mass functions. However, the BH mass function 
has been directly measured only locally, so this comparison 
can be performed only at z = 0, and can not be used to glean 
information on the evolution of these two parameters. 

The left panel of Figure compares the BH mass func- 
tion predicted from the combination of &((j,z) and the mildly 
evolving best-fit Mbh-c relation with a normalization oc (1 + 
z)° 26 (shown as thick curves, with their uncertainty shown in 
gray), and the mass function predicted from the AGN lumi- 



nosity function of SWM assuming a mean A = 0.6 (shown as 
thin curves). 6 

Figure shows that up to z < 3, a constant (non-evolving) 
mean Eddington ratio of A = 0.6 provides a good match be- 
tween the shapes of the accreted BH mass function and the 
one computed from the VDF. At the low-mass end, the VDF- 
based BH mass function starts being incomplete, while at the 
high mass end, a higher intrinsic scatter in the Mbh-c rela- 
tion and/or a more complicated Eddington ratio distribution 
may improve the match. Fully matching the two BH mass 
functions is beyond the scope of this paper (see SWM for fur- 
ther analysis). Our aim here is merely to demonstrate that 
our simple approach also provides hints on the mean Edding- 
ton ratio and its redshift evolution. Similar results are found 
switching to the Hopkins et al. (2007) luminosity function. 
The right panel of Figure0shows that a good match between 
the BH mass functions is recovered on adopting a constant 
A = 1.0. Although systematic uncertainties in the bolomet- 
ric AGN luminosity function preclude tighter constraints on 
the mean Eddington ratio (see SWM for further discussions 
on these issues), it is remarkable that simple models with 
0-5 < A < 1 .0 constant with redshift, can provide a reasonable 
match with the VDF-based BH mass functions. An indepen- 
dent way to constrain the Eddington ratio distribution and its 
evolution with redshift can be derived by matching the halo 
clustering implied by the redshift dependent model BH mass 
function and the observed AGN clustering (Shankar et al., in 
preparation). We have also checked that the same values of A 
provide a good match even at z > 3, however the large uncer- 
tainties associated to the VDF at these high redshifts prevent 
any firm conclusion. 

A somewhat different version of the above exercise was 
performed by Haiman et al. (2007). Under the assumption 
that the duty cycle of quasar activity is short, Haiman et al. 
(2007) matched the instantaneous quasar LF at each redshift 
to the LF predicted from Pvdf(z), plus an assumed constant 
(non-evolving) duty-cycle and Eddington-ratio distribution. 
This approach neglects the BH mass accreted during the lu- 

6 Note that we assumed an initial duty cycle of 0.5 at z = 6; however the 
BH mass function at z < 3.5 becomes independent of this assumption. See 
SWM for further details. 
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FIG. 3. — Confidence levels of 90% and 99% computed assuming x 2 = 
XZm, + 2.30 and \ 2 = X^in + 9.21, respectively, for two parameters in the 
mo del , the radiative efficiency e and the exponent a, where the normalization 
of the Mbh-o" relation evolves as <x (1 +z) a ■ The cross marks the best-fit 
values of e = 0.08 and a = 0.26 corresponding to the minimum Xmm ~ ^ 
(for 8 degrees of freedom). It is evident that a strong redshift evolution in 
the Mbh-o" relation is ruled out at a high confidence level, if the radiative 
efficiency is constant in time. 

minous quasar phases (or at least any corresponding varia- 
tion of the "quasar light-curve" caused by the growth in BH 
mass), and places a constraint directly on the relation between 
quasar luminosity L and host velocity dispersion a. While the 
L — a relation is essentially a convolution of the Eddington ra- 
tio distribution with the Mbh-c relation, this approach cannot 
be used to study these two relations separately. Nevertheless, 
Haiman et al. (2007) found no evidence for any evolution in 
the L — a relation with redshift; their fits to the quasar LF are 
consistent with a constant 0.3 <> A < 0.5 combined with a non- 
evolving Mbh-c relation. Since the Mbh-c relation is indeed 
found here, independently, to be non-evolving, this breaks the 
degeneracy in the result of Haiman et al. (2007) and also re- 
quires that the evolution in the Eddington ratio distribution be 
modest. 
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FIG. 4. — Same format as Figure [2] with parameters a = 1 and e ~ 
0.08.0.10 for the SWM and Hopkins et al. (2007) luminosity functions, re- 
spectively. 

4. DISCUSSIONS 

4. 1 . VARYING THE MODEL ASSUMPTIONS 

We have studied more complicated scenarios where we also 
allow for the scatter and/or the slope of the Mbh-c relation to 
increase with redshift. For example, steadily increasing the 



slope from 3.83 to, say, 5.5, at fixed scatter of 0.22 dex, still 
implies a ~ 0.30. The left panel of Figure [6] shows instead 
the comparison between the p,(z) and the Pvdf(z), assuming 
the scatter increases with redshift from -q = 0.22 at z = to 
?y = 0.40 at z = 5.7, the highest redshifts probed by our sam- 
ple. It can be seen that the best-fit model requires a = 0.15, 
even lower than what reported in Figure|2] This is expected as 
these models tend to increase the BH mass density associated 
with the VDF at a given redshift, implying an even milder de- 
gree of evolution in the Mbh-c normalization. The right panel 
of Figure [6] compares the implied mass functions predicted 
by the same ?7(z)-model and by the AGN luminosity func- 
tion. Similarly to the best-fit model discussed in § [3] a good 
match can be recovered assuming a constant A = 0.6. More- 
over, steadily increasing the intrinsic scatter from 0.22 to 0.4 
dex significantly improves the match between the VDF- and 
AGN-based BH mass functions at both the high and low-mass 
ends. 

We have also explored different models for the redshift evo- 
lution in the Mbh-c relation. For example, a model in which 
a quadratic term of the type 5 x log(l +z) 2 is added in equa- 
tion (f2]l (see Wyithe 2004 for a similar test applied to the local 
A^bh-c relation), produces a good match between Pvdf(z) and 
p, (z) if a ~ S ~ 0.1, with \ 2 ~ 2.5, comparable to our best-fit 
model. We therefore conclude that although the choice for the 
redshift evolution model in equation|2]is not unique, alterna- 
tive solutions will still provide similar constraints on the net 
amount of allowed evolution. 

4.2. THE SCATTER IN THE M BH - <r RELATION 

At each redshift, our approach assumes that new BHs are 
formed with a mass that is tightly imposed by the velocity 
dispersion of the host galaxy. Increasing the normalization 
and/or scatter of the Mbh-c relation at high redshift therefore 
induces in the local universe a finite spread in BH mass at 
fixed velocity dispersion. If the evolution is steep, this can 
exceed the observed scatter r/ < 0.22 dex. Each panel in Fig- 
ure |7]plots as solid lines the median BH mass distribution of 
100 Monte Carlo simulations corresponding to a given bin of 
velocity dispersion, as labeled. The BH masses are derived 
from the redshift-dependent Mbh-c relation where the red- 
shifts are randomly extracted from the age distribution Pij(a) 
competing to each velocity dispersion. The long-dashed lines 
refer to the Gaussian distribution with r\ = 0.22 measured by 
Tundo et al. (2007). Both our models with constant (upper 
panels) or evolving (lower panels) scatter still produce at z = 
BH mass distributions at fixed velocity dispersion comparable 
to what is observed. The small off-set in the Gaussian dis- 
tributions predicted by our simulations with respect to those 
observed is induced by the sampling of higher redshift, more 
massive BHs. Models characterized by stronger redshift evo- 
lution with a > 0.3 will then evidently predict a scatter in the 
local Mbh-c relation much larger than what is actually ob- 
served. Mergers are then required to be a significant compo- 
nent in the evolution of the BH population in these models, as 
the Monte Carlo simulations performed by Peng (2007) show 
that random BH mergers will tighten the relations between 
BH and host galaxy masses at late times. However, frequent 
mergers may, on the other hand, predict too many massive 
BHs with respect to those seen in the local universe (see Fig- 
ure 13 in SWM). 

4.3. SYSTEMATIC UNCERTAINTIES IN THE METHOD 
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The main result of this paper is shown in the right panel 
of Figure [2] which demonstrates that a good match between 
Pvdf(z) and p,(z) can be achieved by assuming a mild red- 
shift evolution in the Mbh-c relation with a < 0.3. These 
results are based on the age distributions pjj derived from 
Bernardi et al. (2006). However, MOPED-based galaxy ages 
are, on average, larger at fixed velocity dispersion, predicting 
a flatter dependence $(c,z) as a function of redshift z (see 
Figure Q]). This will correspondingly flatten Pvdf(z) vs red- 
shift, and decrease the best-fit a, therefore requiring an even 
milder redshift evolution in the Mbh-c relation. A null evolu- 
tion in the Mbh-c relation is expected in basic AGN feedback 
models (e.g., Silk & Rees 1998), in which a tight correlation 
derives by imposing equilibrium between the energy released 
by the central BH, and the gas binding energy, linked to the 
velocity dispersion. 

Bernardi et al. (2007), Graham (2007), and Shankar & Fer- 
rarese (2008) have discussed selection biases in the available 
sample of BHs that may induce systematic uncertainties in the 
determination of the local BH mass function. However, our 
conclusions are not affected by these uncertainties, because 
a change in the local BH mass density would be absorbed in 
the radiative efficiency e (i.e. e would be modified, to match 
Pvdf(z) and p,(z) at z = 0, but a would not change). By the 
same token, our results are only weakly dependent on whether 
or not the bulges of spirals are included in the estimate of 
the local BH mass function (a weak dependence arises only 
because the addition of the spiral bulges slightly skews the 
age-distribution of the total population to younger ages; this 
becomes increasingly less important toward higher redshifts, 
where a progressively smaller fraction of the total BH mass 
density is contributed by the low-er galaxies). 

Likewise, uncertainties in redshift-independent bolometric 
corrections do not alter our conclusions. The bolometric cor- 
rection adopted in SWM is lower by ~ 30% with respect to 
the one used by Hopkins et al. (2007), but the sole effect 
of this difference is to yield a proportionally smaller value of 
the mean radiative efficiency to recover the match between 
Pvdf(z) and p,(z) at z = (see left panel of Figure©. More- 
over, the break luminosity and bright-end slopes of the SWM 
and Hopkins et al. (2007) luminosity functions are somewhat 
different (see Figure 4 in SWM). Nevertheless, within uncer- 
tainties, the resulting BH accretion histories obtained from the 
two luminosity functions have a similar behavior with red- 
shift, both placing the same constraint a < 0.3 for the evolu- 
tion in the normalization of the Mbh-c relation. 

On similar grounds, if we assume that the BHs in our 
sample radiate at even lower luminosities than the L m ; n con- 
sidered in the integral of equation (0, our results do not 
change. For example, lowering the minimum luminosity to 
logL m i n /ergs _1 = 42, the cumulative emissivity of AGNs in- 
creases by about ~ 30% at all redshifts yielding a very sim- 
ilar behavior with time. Therefore, a proportionally higher 
radiative efficiency plugged into equation (01 keeps the good 
match with the Pvdf(z) and specifically we find that the value 
of a — 0.16 yields a^ 2 ~ 2.4 for 8 degrees of freedom. 

Our conclusions about the (lack of) evolution in the normal- 
ization of the Mbh-c relation, in general, are more dependent 
on redshift-dependent effects. For example, if the bolomet- 
ric correction increased to high z (or, e.g., if obscuration were 
more significant at higher redshift), this would again further 
decrease our favored mild positive redshift evolution in the 
Mbh-c normalization. Likewise, evolution in the mean radia- 



tive efficiency and/or the assumed scatter in the Mbh-c rela- 
tion would modify our results, in the sense that our predicted 
evolution would be milder if either increased toward high z. 
In principle, to allow for a stronger evolution in the Mbh-c 
relation the radiative efficiency must significantly decrease at 
z > 3 to boost the accreted mass density at fixed AGN lumi- 
nous density. However, we have checked that e must then 
rapidly increase at lower redshifts in order not to overpro- 
duce the local BH mass density. More quantitatively, if we set 
e - 0.05 at z > 3, then it must be that e > 0.05 x [7/(1 + z)] a5 
at lower redshifts. Such an evolution in e is not enough to 
allow for a strong variation in the Mbh-c relation. We found 
that Pvdf(z) can match the p. (z) implied by the e(z)-model 
if a ~ 0.3, which is close to our best-fit model. On other 
grounds, as recently shown by Shankar et al. (2008b), a too 
low radiative efficiency at high redshifts seems to be disfa- 
vored by BH accretion models which simultaneously repro- 
duce the strong quasar clustering measured at z — 3 — 4 in 
SDSS by Shen et al. (2007a), the mean Eddington ratio of 
A > 0.5, measured by Shen et al. (2007b) for the same quasar 
sample, and the high redshift quasar luminosity function (e.g., 
Richards et al. 2006; Fontanot et al. 2007; Shankar & Mathur 
2007). 

4.4. EVOLVING THE MAGORRIAN RELATION 

Most of the results from other groups discussed in §[T]focus 
on the ratio between black hole mass and stellar mass. The lat- 
ter may settle on longer timescales with respect to the galaxy 
velocity dispersion, the amplitude of which is linked to the 
central potential well which grows faster than the overall evo- 
lution of the halo (Zhao et al. 2003). In order to get some hints 
on the actual evolution of the Mbh-Mstar relation with red- 
shift, we have converted the galaxy stellar mass function into 
a BH mass function assuming the Mbh-Mstar ratio evolv- 
ing as (1 +z)», We have used the recent near-infrared stel- 
lar mass function by Perez-Gonzalez et al. (2008), well con- 
strained within < z < 3 and 10 < log M*/M < 12. We have 
then converted the latter into a BH mass function by assum- 
ing that, on average, about 0.7 x 10~ 3 (1 + z)^ (e.g., Magor- 
rian et al. 1998; Marconi & Hunt 2003) of the total stellar 
mass is locked up in spheroids and is associated to the central 
black hole, with a Gaussian scatter around the mean of 0.3 dex 
(e.g., Hiiring & Rix 2004). In this case, we find that ( < 0.3 
is a necessary condition for the BH mass density to be con- 
sistent with the accreted mass from AGNs, the latter derived 
assuming a fixed value of the radiative efficiency. This re- 
sult is in agreement with the degree of evolution discussed in 
§[3]found by evolving the Mbh-c relation. Although these re- 
sults are in reasonable agreement with other works (Marconi 
et al. 2004; Merloni et al. 2005; De Zotti et al. 2006; SWM; 
Merloni & Heinz 2008), uncertainties on the lower limit of 
the stellar mass function and/or on the true fraction of stellar 
mass associated to BH growth at any time, make this method 
less reliable than the one based on velocity dispersion, and we 
therefore do not pursue it further. 

4.5 . THE IMPACT OF MERGERS 

So far we have neglected mergers in our calculations. Ma- 
jor mergers between massive galaxies do occur, although re- 
cent work has suggested the galaxy merger rate may be lower 
than previously thought. Drory & Alvarez (2008) compared 
the time variation in the stellar mass function with the evolu- 
tion implied by the star formation rate alone, concluding that 
galaxies with stellar masses above 10 u Mq undergo at most 
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FIG. 5. — Comparison between the mass functions predicted from <E>(o\z) convolved with the Mbh-o" relation whose normalization evolves as oc (1 +z) ' 26 
(thick curves, with the uncertainty shown shaded in grey), and the mass function predicted from the AGN luminosity functions of Shankar et al. (2008; left panel) 
and Hopkins et al. (2007; right panel). A constant Eddington ratio of A = 0.6, 1.0 has been assumed for computing the accreted mass functions in the left and 
right panels, respectively. It can be seen that the choice of a single constant Eddington ratio provides a good match to the velocity dispersion-based black hole 
mass functions, at least around the peak of the distributions. A good match is also found extending the comparison up to z > 3, however the large uncertainties 
at these redshifts prevent any firm conclusion. 



one major merger since z ~ 1.5, in agreement with the results 
of Bell et al. (2007). Lotz et al. (2006) find evidence for an 
even lower merger rate since z ~ 1 from the DEEP2 survey. 
Most importantly, however, a significant rate of major merg- 
ers would strengthen our conclusions. In velocity dispersion 
space, collisionless major mergers do not significantly affect 
the final a. For example, in a dry merger of comparable-mass 
galaxies with mass Mi and M 2 and corresponding velocity 
dispersions <j\ and 02, the resulting galaxy will have a veloc- 
ity dispersion cr 2 ~ \M\o\ +M 2 aj]/(M l +M 2 ) <max(er 2 ,o- 2 ) 
(e.g., Ciotti et al. 2007). Therefore, if the masses of the two 
galaxies are comparable, the final a will be close to the veloc- 
ity dispersion of the progenitors. 

Dry mergers would then double the number of galaxies we 
predict at fixed a. Every dry major merger would in fact 
"split" the galaxy into two (going back in time), adding an 
extra galaxy with the same a, compared to our present pre- 
dictions (which neglect mergers). In turn, this implies that 
the merger-free VDF $(cr,z) computed above for z > is a 
lower limit to the true abundance of galaxies with velocity 
dispersion a at redshift z. The associated BH mass density 
p» (z) will consequently also be underestimated at redshift z. 
The inclusion of any mergers then predicts a larger BH mass 
density at fixed a; to compensate for this increase, a lower 
normalization of the Mbh-c relation at z = 1 — 2 is then re- 
quired, which strengthens our conclusions that large a values 
are excluded by the match between Pvdf(z) and p,(z) (see 
the left panel in Figure |2). 

However, if dissipation played a non-negligible role during 
the evolution of the galaxy (either as a result of mergers, or in 
isolation), then the velocity dispersion may increase with time 
from the epoch of first collapse. To mimic such effects, we al- 
low all velocity dispersions to increase at higher redshifts as 
<j(z) = u(0) x (1 +z)~ 7 . Most probably this evolution is mass 
and/or velocity-dependent, nevertheless this approach will be 
able to set interesting constraints on the mean variation of a. 
Also, any estimate for 7 should here be considered as a lower 
limit to the actual evolution of a, as we neglect the still poorly 
understood increase in galaxy number density due to possi- 
ble galaxy mergers. In Figure [8] we show the main results 
for a model with 7 = 0.23 and a — 1.5. Note that with the 
adopted scaling Mbh oc a" with j3 ~ 4, we expect, a degener- 
acy between 7 and a given approximately by7«a//3«a/4 



(although the degeneracy is modified slightly by the assumed 
scatter and age-spread of BHs at a given z and a). We find 
that the downsizing evolution in this case is canceled out (up- 
per left panel), as all galaxies are now pushed to lower and 
lower a at higher redshifts. 

The evolving number density in this model, shown in the 
upper right panel of Figure [8] seems to be at variance with 
the number density evolution of early-type galaxies within 
< Z < 1 inferred from DEEP2 by Faber et al. (2007; solid 
points in the same Figure). However, lacking a clear under- 
standing of how the &((T,z) should evolve in the presence of 
mergers, this model cannot be ruled out, although some incon- 
sistencies can already be pointed out. From the Figure it can 
be seen that the number density of the massive red galaxies in 
DEEP2 at the intermediated redshifts of 0.5 < z < 1 , should be 
matched with galaxies characterized by a velocity dispersion 
a < 100 kms -1 . Alternatively, the DEEP2 number densities 
could be matched with the number density of galaxies with 
larger velocity dispersions if mergers were a significant com- 
ponent in the evolution of these galaxies, thus significantly 
increasing their number density at higher redshifts. However, 
the latter hypothesis may contradict independent results (e.g., 
Lotz et al. 2006). The dissipative model described here also 
predicts a mean velocity dispersion about flat out to z ~ 2 and 
slightly decreasing at higher redshifts, as shown in the lower 
left panel of Figure [8] 

The main achievement of this model is the good match 
between Pvdf(z) and p»(z) even if a strong evolution in the 
Mbh-c relation has been assumed (a = 1.5), as shown in the 
lower right panel of Figure|8] This model is characterized by a 
significant dissipative phase in the evolution of typical early- 
type galaxies, which could represent an interesting constraint 
for galaxy evolution models and it can in principle be tested 
through hydrodynamical simulations, which we plan to do in 
future work. 

On the other hand, a major problem with the dissipative 
model is represented by its implied Eddington ratio distribu- 
tion. We in fact find that the strong increase in the Mbh-c nor- 
malization at higher redshifts requires a significant decrease, 
by up to a factor of a few, in the mean Eddington ratio A(z) to 
keep the match between the BH mass functions at z > 2 shown 
in Figure|5] The latter behavior of A(z) is at variance with sev- 
eral works which actually claim an almost constant or proba- 
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FIG. 6. — Left panel: comparison between the p. (z) and the Pvdf(z) computed for a model in which the $(ct,z) is convolved with aM^-a relation the scatter 
of which increases with redshift from r] = 0.22 to 0.40 . Right panel: comparison between the mass functions predicted from the same model and from the AGN 
luminosity functions with A = 0.6. 
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FIG. 7. — In each panel the solid histograms show the mean of a set of Monte Carlo simulations which compute from our models the expected distribution of 
black hole masses at fixed velocity dispersion, as labeled (see text for details). The long-dashed lines refer to the Gaussian distribution with r\ = 0.22 measured 
by Tundo et al. (2007). Both our models with constant (upper panels) or evolving (lower panels) scatter still produce at z = black hole mass distributions at 
fixed velocity dispersion comparable to what is observed. Models characterized by stronger redshift evolution will then evidently predict a scatter in the local 
Mftn-a relation much larger than what is observed. 



bly increasing A(z) at higher redshifts (e.g., McLure & Dun- 
lop 2004; Shankar et al. 2004; Vestergaard 2004; Kollmeier 
et al. 2006; Netzer & Trakhtenbrot 2007; SWM; Shen et al. 
2008b; Shankar et al. 2008b). 

4.6. COMPARISON WITH PREVIOUS WORKS 

The relatively mild Mbh-c redshift evolution inferred from 
our approach may seem in apparent disagreement with some 
recent independent studies. As briefly mentioned in §Q] Treu 
et al. (2007) and Woo et al. (2008) have randomly compiled 
from the SDSS Data Release 4 a sample of about 20 Seyferts 
galaxies in the redshift range 0.37 < z < 0.57. Their results, 
shown as open circles in Figure [9] are compared with those 
of Shen et al. (2008, shown as filled circles), who estimated 
the Mbh-c relation for a larger sample of active galaxies up to 
z = 0.452. While the latter claim that no significant evolution 
in the Mbh-c relation is detectable from their sample, Woo 
et al. (2008) confirm the results by Treu et al. (2007) that a 



significant increase of ~ 0.2 dex in BH mass at fixed velocity 
dispersion must occur within z = and z ~ 0.5. Our best-fit 
model, shown at redshifts z = and z = 0.5 with long-dashed 
and solid lines respectively, shows no strong evolution within 
this redshift range and it is in reasonable agreement with both 
samples. A significant discrepancy is noticeable with respect 
to the Woo et al. (2008) results for velocity dispersions log 
<r/kms~ l < 2.3. However, systematic uncertainties may af- 
fect these estimates; for example as also discussed by Woo 
et al. (2008), especially in galaxies with lower BH mass, the 
host galaxy contribution to the 5 100 A luminosity may lead 
to an overestimation of the true BH mass. Overall, given the 
systematics and biases which affect these kind of studies (e.g., 
Lauer et al. 2007), we do not find strong evidence for a dis- 
agreement between these works and our results. For the same 
reasons, we do not attempt comparisons with the results ob- 
tained from higher redshift studies. 
Merloni et al. (2004) compared the accreted BH mass den- 
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FIG. 8. — Results for a model in which we allow the mean galaxy velocity dispersion to decrease at higher redshifts to mimic the effects of prolonged wet 
activity in the host galaxies since their formation epoch. We assume an evolution of the type a(z) = c(0) X (1 + z) * with 7 = 0.23. The downsizing effect is 
canceled out (upper left) and the match with DEEP2 number density evolution is significantly worsen (upper right). At variance with our previous results this 
model predicts a mean velocity dispersion about flat out to z ~ 2 (lower left) and, more important, a good match between pvdf(z) and p, (z) with a = 1 .5 (lower 
right). 



sity in AGNs with the cosmological global star formation rate 
density (see also Haiman et al. 2004). Although their conclu- 
sions depend on additional assumptions about the fraction of 
the star forming galaxies which are linked to BH growth at a 
given redshift, irrespective of the adopted value of the radia- 
tive efficiency their best-fit relation yields a « 0.5, somewhat 
higher, but still consistent, with the value found here, and they 
also rule out a > 1.2 at a high confidence level. Hopkins et 
al. (2006) also describe a model-independent integral con- 
straint that defines an upper limit to the allowed degree of 
evolution in the ratio of BH mass to host galaxy luminosity 
or mass, as a function of redshift. By comparing the AGN 
density with the luminosity and mass functions in different 
bands from redshifts z = — 2, they rule out at > 6a a BH- 
host galaxy mass ratio significantly larger at high redshifts 
than locally. Cattaneo & Bernardi (2003) combined a rela- 
tion between mean age and velocity dispersion, derived from 
a sample of SDSS local early-type galaxies, with the Sheth et 
al. (2003) local VDF. By assuming a redshift independent 
mean Eddington ratio, radiative efficiency and obscuration 
correction, they were then able to reproduce the AGN optical 
and X-ray luminosity functions. As mentioned above, similar 
calculations have been performed recently by Haiman et al. 
(2007), whose results imply, assuming a non-evolving Mbh- 
g relation, that in order to reproduce the bolometric quasar 
luminosity function, the quasars must shine at a mean sub- 
Eddington regime of A = 0.5 that is approximately constant 
with time. This conclusion was confirmed by the indepen- 



dent estimates of SWM. However, the works by Cattaneo & 
Bernardi (2003) and Haiman et al. (2007) can only constrain 
the combination of the Eddington ratio distribution and the 
A^bh-c relation, while our approach here can simultaneously 
constrain the mean accretion histories of BHs and their host 
galaxies, and the mean Eddington ratio of BHs at all times. A 
further difference is that the analysis of Haiman et al. (2007) 
can constrain the quasar lifetime, while the results here rely on 
the comparison between time-integrated quantities, and are 
strictly independent of the quasar lifetime. 

5. CONCLUSIONS 

In this work we combined the local VDF with the stellar 
age distributions estimated by Bernardi et al. (2006), to com- 
pute the VDF at higher redshifts, <&(ct,z). In agreement with 
previous work, we find statistical evidence for downsizing, 
whereby the stellar populations in galaxies with higher ve- 
locity dispersions formed earlier, irrespective of the specific 
model we adopt for computing the galactic ages. We then 
computed the BH mass function associated with $(ct,z) at 
each redshift z, through a BH mass - velocity dispersion rela- 
tion whose normalization was allowed to evolve with redshift 
as oc (1 + z) a ■ Our main underlying assumptions are that most 
of the growth of the central BH occurs simultaneously (within 
± 1 Gyr) with the formation of the host's potential well, and 
that the measured stellar ages represent this formation time 
to within a similar accuracy. The BH mass density Pvdf(z) 
inferred from the VDF can then be compared with the accu- 
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FIG. 9. — Our best-fit model for the M^H-cr relation is plotted at z = and 
z = 0.5, as labeled, and compared with recent data from Treu et al. (2007) 
and Woo et al. (2008), shown with open symbols, and Shen et al. (2008), 
shown with filled symbols. 

mulated BH mass density implied by the time-integral of the 



AGN luminosity function, p. (z) . We find significant evidence 
that the match between pvdf(z) andp,(z) implies a relatively 
mild redshift evolution, with a < 0.30, and with values of 
a > 1.3 excluded at 99% confidence. If a positive redshift 
evolution stronger than a > 1 were to be confirmed indepen- 
dently in the future, then this would be a robust indication that 
dissipative processes played a significant role in galaxy evo- 
lution, resulting in an increase in the velocity dispersion of 
the spheroid components of individual galaxies with cosmic 
time. However, we also find evidence that a dissipative model 
predicts a mean Eddington ratio decreasing with increasing 
redshift, at variance with several independent studies. 
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